---
title: "R21.Expropriation_Robust"
---

# Load Packages
```{r}
library(tidyverse)
library(estimatr)
```

# Load Data
```{r}
load("Data/data_final.Rdata")
```

# Plot Relationship between Expropriatio and Communist Vote Share in 1975 and Distance to Nearest PCP Office (App. Figures 8aii)
```{r}
data_final %>% 
  ggplot(aes(x=pr_ha_ord, y=Communist_1975)) + 
  geom_point(alpha = 0.25) +
  geom_smooth(col = "black") + 
  geom_smooth() + 
  ylab("Communist Vote Share, 1975") + xlab("% Hectares Expropriated") +
  theme_light() 

ggsave(plot = last_plot(),
             "Figures/Afg8a2.1_TrExpPCP75.pdf")

data_final %>% 
  ggplot(aes(x=pr_ha_ord, y=MinDistPCP_Map/1000)) + 
  geom_point(alpha = 0.25) +
  geom_smooth(col = "black") + 
  geom_smooth() +
  ylab("Distance to Nearest PCP Worker Office") + xlab("% Hectares Expropriated") +
  theme_light() 

ggsave(plot = last_plot(),
             "Figures/Afg8a2.2_TrExpPCPOffice.pdf")
```

# Model Relationship between Expropriation and Communist Vote Share in 1975 and Distance to Nearest PCP Office (App. Table 8aii)
```{r}
models<-list()

models[[1]]<-lm_robust(pr_ha_ord ~ Communist_1975 + 
                         ln_Confs_All + as.factor(DistN) +ln_PopDens_70 +Cities_10km +PrRecH_Adm3 +Isolados1950 +PrWheatAr_1972 +LandIneqRat +TerrainSlopeIndex +MedAltitude, 
                       data_final) %>% summary()
models[[2]]<-lm_robust(pr_ha_ord ~ MinDistPCP_Map + 
                         ln_Confs_All + as.factor(DistN) +ln_PopDens_70 +Cities_10km +PrRecH_Adm3 +Isolados1950 +PrWheatAr_1972 +LandIneqRat +TerrainSlopeIndex +MedAltitude, 
                       data_final) %>% summary()

output_covs <- data.frame(seq = 1:1, 
                          variable = c("Communist Vote Share, 1975", "Distance to Nearest PCP Office"), 
                          outcome = c("Expropriation", "Expropriation"),
                          coef = NA, se = NA, n = NA)

for (i in 1:length(models)){
  output_covs[i, "coef"]<-models[[i]]$coefficients[2, 1]
  output_covs[i, "se"]<-models[[i]]$coefficients[2, 2]
  output_covs[i, "n"]<-models[[i]]$df.residual
}

ggplot(output_covs, aes(x = variable)) +
  geom_hline(aes(yintercept = 0), col = "red") +
  geom_point(aes(y = coef), position = position_dodge(width = 1)) + 
  geom_linerange(aes(ymin=coef - (qt(0.975, n)*se), ymax=coef + (qt(0.975, n)*se)), position = position_dodge(width = 1), size = 0.5) +
  geom_linerange(aes(ymin=coef - (qt(0.95, n)*se), ymax=coef + (qt(0.95, n)*se)), position = position_dodge(width = 1), size = 0.75) +
  coord_flip() +
  ylab("Effect of IVs on % Hectares Expropriated") + xlab("Independent Variables") +
  theme_light() 

ggsave(plot = last_plot(), 
              "Figures/Afg8a2.3_BalExprop.pdf", 
       height = 2, width = 7)

```